#ifndef ARRAY_HH
#define ARRAY_HH

#include <algorithm>

//////////////////////////////////
//	constructors
//////////////////////////////////

template <typename T, unsigned int DIM>
array<T, DIM>::array( const array<T, DIM>& ref_arr )
{
	*this = ref_arr;
}

template <typename T, unsigned int DIM >
array<T, DIM>& array<T, DIM>::operator=( const array<T, DIM>& ref_arr )
{
	int num_elem = (*ref_arr.data_p_).size();
	data_p_ = std::shared_ptr<storage_type>( new storage_type( num_elem ) );
	std::copy( (*ref_arr.data_p_).begin(), (*ref_arr.data_p_).end(), (*data_p_).begin() );
	
	storage_ = ref_arr.storage_;
	shape_ = ref_arr.shape_;
	stride_ = ref_arr.stride_;
	p_zero_offset_ = ref_arr.p_zero_offset_;
	return *this;
}

template <typename T, unsigned int DIM>
array<T, DIM>::array( array<T, DIM>&& arr )
{
	*this = std::move( arr );
}
	
template <typename T, unsigned int DIM>
array<T, DIM>& array<T, DIM>::operator=( array<T, DIM>&& ref_arr )
{
	data_p_ = std::move( ref_arr.data_p_ );
	storage_ = std::move( ref_arr.storage_ );
	shape_ = std::move( ref_arr.shape_ );
	stride_ = std::move( ref_arr.stride_ );
	p_zero_offset_ = std::move( ref_arr.p_zero_offset_ );
	return *this;
}

template <typename T, unsigned int DIM>
array<T, DIM>::array( std::initializer_list<int> shape )
	: storage_( array_storage<DIM>::c_array() ), stride_(1)
{
	SPX_ASSERT( shape.size() == DIM || shape.size() == 1, "length of dimensions does not match" );
	
	shape_ = (shape.size() == 1) ? 
		tiny_vector<int, DIM>( *(shape.begin()) ) : tiny_vector<int, DIM>( shape );

	create_new_array();
	array_helper::init_pzero_offset( *this, p_zero_offset_ );
}

template <typename T, unsigned int DIM>
array<T, DIM>::array( tiny_vector<int, DIM> shape, array_storage<DIM> storage)
	: storage_(storage), shape_(shape), stride_(1)
{
	SPX_ASSERT( shape.size() == DIM, "length of dimensions does not match" );
	
	create_new_array();
	array_helper::init_pzero_offset( *this, p_zero_offset_ );
}

template <typename T, unsigned int DIM>
array<T, DIM>::array( tiny_vector<range, DIM> range, array_storage<DIM> storage )
	: storage_(storage), stride_(1)
{
	tiny_vector<int, DIM> base;
	for( unsigned int dim = 0; dim < DIM; dim++ )
	{
		base[dim] = range[dim].start();
		shape_[dim] = range[dim].length();
	}
	storage_.set_base( base );
	
	create_new_array();	
	array_helper::init_pzero_offset( *this, p_zero_offset_ );
}

//! for expression
template <typename T, unsigned int DIM>
template <typename expr_type>
array<T, DIM>::array( expression<expr_type> expr )
{
    // Determine extent of the array expression
	tiny_vector<int, DIM> lbound, shape, ordering;
	tiny_vector<bool, DIM> in_ordering = false;
	
	int j = 0;
	for( int i = 0; i < DIM; ++i)
	{
		lbound(i) = expr.lbound(i);
		shape(i) = expr.shape(i);
		auto orderingj = expr.ordering(i);
		static_assert( !std::is_same<decltype(orderingj), unbound>::value, "unbound array expression" );
		// unique value in ordering array
		if( orderingj < DIM && !in_ordering( orderingj ) )
		{ 
			in_ordering( orderingj ) = true;
			ordering( j++ ) = orderingj;
		}
	}

	// It is possible that ordering is not a permutation of 0,...,N_rank-1.
	// In that case j will be less than N_rank. We fill in ordering with the
	// usused values in decreasing order.
	for( int i = DIM - 1; j < DIM; ++j ) 
	{
		while( in_ordering(i) )
			--i;
		ordering(j) = i--;
	}
		
	array<T, DIM> A( shape, array_storage<DIM>( ordering, lbound ) );
	*this = std::move(A);
	*this = expression<decltype(expr)>(expr);
}

//! for slicing
template <typename T, unsigned int DIM>
struct slice_callback
{
	slice_callback( array<T, DIM>& ref_arr ) : ref_arr_( ref_arr ) {}
	void operator()( int dim, range r ) 	{ ref_arr_.slice( dim, r ); }
	void operator()( int dim, range_all r ) { ref_arr_.slice( dim, range( ref_arr_.lbound(dim), ref_arr_.ubound(dim) ) ); }
	void operator()( int dim, int r ) 		{ ref_arr_.slice( dim, range( r, r ) ); }
	array<T, DIM>& ref_arr_;
};

template <typename T, unsigned int DIM>
template <typename... range_types>
array<T, DIM>::array( array& ref_arr, std::tuple<range_types...> new_slice )
{
	reference( ref_arr );
	static_assert( array_helper::range_types_checker<array, range_types...>::valid, "invalid argument");
	slice_callback<T, DIM> f( *this );
	tuple_for_each( f, std::forward<std::tuple<range_types...>>(new_slice) );
}

//! for exisiting memory block
template <typename T, unsigned int DIM>
array<T, DIM>::array( std::shared_ptr<storage_type> data, tiny_vector<int, DIM> shape, array_storage<DIM> storage )
	: storage_(storage), shape_(shape)
{
	stride_ = storage_.compute_storage_stride( shape_ );	
	array_helper::init_pzero_offset( *this, p_zero_offset_ );
	data_p_ = data;
}

template <typename T, unsigned int DIM>
array<T, DIM>::array( std::shared_ptr<storage_type> data, tiny_vector<int, DIM> shape, 
	tiny_vector<int, DIM> stride, int p_zero_offset, array_storage<DIM> storage )
	: data_p_(data), storage_(storage), shape_(shape), stride_(stride), p_zero_offset_(p_zero_offset)
{ }

//////////////////////////////////
//	multi-components
//////////////////////////////////
template <typename T, unsigned int DIM>
template <template<class multi_comp_T> class multi_component_checker>
array<typename multi_component_checker<T>::component_type, DIM>
array<T, DIM>::extract_component( unsigned int component_id )
{
	typedef typename multi_component_checker<T>::component_type comp_T;
	typedef array<comp_T, DIM> return_arr_T;
	typedef typename return_arr_T::storage_type new_storage_T;
	static constexpr unsigned int num_comps = multi_component_checker<T>::num_components;

	SPX_PRE_CONDITION( component_id < num_comps );
	
	new_storage_T* re_p = reinterpret_cast<new_storage_T*>( data_p_.get() );
	std::shared_ptr<new_storage_T> sh_pt( data_p_, re_p );
	tiny_vector<int, DIM> new_stride = stride_;
	for( unsigned int dim = 0; dim < DIM; dim++ ) 
		new_stride(dim) *= (int)num_comps;
	int new_p_zero = p_zero_offset_ + component_id;
	return_arr_T new_arr{ sh_pt, shape_, new_stride, new_p_zero, storage_ };
	
	return std::move(new_arr);
}

template <typename T, unsigned int DIM>
array<typename multi_component_traits<T>::component_type, DIM>
array<T, DIM>::operator []( unsigned int component_id )
{
	return extract_component( component_id );
}

//////////////////////////////////
//	methods
//////////////////////////////////

template <typename T, unsigned int DIM>
void array<T, DIM>::create_new_array()
{
	int total_size = 1;
	for( auto v : shape_ )	total_size *= v;
	data_p_ = std::shared_ptr<storage_type>( new storage_type( total_size ) );
	stride_ = storage_.compute_storage_stride( shape_ );
	p_zero_offset_ = 0;
}

template <typename T, unsigned int DIM>
void array<T, DIM>::reference( array<T, DIM>& ref_arr )
{
	data_p_ = ref_arr.data_p_;
	storage_ = ref_arr.storage_;
	shape_ = ref_arr.shape_;
	stride_ = ref_arr.stride_;
	p_zero_offset_ = ref_arr.p_zero_offset_;
}

template <typename T, unsigned int DIM>
void array<T, DIM>::slice( unsigned int dim, range r )
{
	auto new_first  = r.start();
	auto new_stride = r.stride();
	auto new_length = r.length();
	
	p_zero_offset_ += ( new_first - lbound(dim) * new_stride) * stride_[dim];
	shape_[dim] = new_length;
	stride_[dim] *= new_stride;	
}

template <typename T, unsigned int DIM>
const T* array<T, DIM>::data_first() const
{
	int p = p_zero_offset_;
	auto id0 = lbound();
	for( unsigned int dim = 0; dim < DIM; dim++ )
		p += id0[dim] * stride_[dim];
	return &(*data_p_.get())[p];
}

template <typename T, unsigned int DIM>
T* array<T, DIM>::data_first()
{
	return const_cast<T*>( 
		static_cast<const array<T, DIM>*>(this)->data_first() );
}

template <typename T, unsigned int DIM>
array<T, DIM>& array<T, DIM>::operator=( std::initializer_list<T> vals)
{
	auto iter1 = begin();
	auto iter2 = vals.begin();
	for( ; iter1 != end() && iter2 != vals.end(); ++iter1, ++iter2 )
		*iter1 = *iter2;
	return *this;
}

///////////////////////////////////////////////
//	subscripting / slicing / resize
///////////////////////////////////////////////

template <typename T, unsigned int DIM>
template <typename... ID_types>
typename array_helper::array_slice<array<T, DIM>, ID_types...>::return_type
array<T, DIM>::operator()( ID_types... id )
{
	return array_helper::array_slice<array_type, ID_types...>::get( *this, id... );
}

template <typename T, unsigned int DIM>
template <typename... ID_types>
typename array_helper::array_slice<const array<T, DIM>, ID_types...>::return_type
array<T, DIM>::operator()( ID_types... id ) const
{
	return array_helper::array_slice<const array_type, ID_types...>::get( *this, id... );
}

template <typename T, unsigned int DIM>
template <typename... ID_types>
void array<T, DIM>::resize( ID_types&&... args )
{
	array_helper::array_resize<array_type, ID_types...>::resize( *this, std::forward<ID_types>(args)... );
	create_new_array();
}

//////////////////////////////////////////////////////
//	operator overloading -- binary assign
//	-- On const value
//////////////////////////////////////////////////////

#define DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( name, op )		\
template <typename T, unsigned int DIM>			\
array<T, DIM>& array<T, DIM>::operator op( const T const_val )	\
{									\
	using binary_assign_ops::name;					\
	iterator_on_one_element<T> iter_wrap( const_val );		\
	/*array_binary_apply<name>::go( begin(), end(), iter_wrap );	*/\
	array_helper::array_evaluator<name>::evaluate( *this, iter_wrap ); \
	return *this;							\
}

DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( equal, = );
DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( plus_equal, += );
DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( minus_equal, -= );
DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( multiply_equal, *= );
DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( divide_equal, /= );
DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( mod_equal, %= );
DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( xor_equal, ^= );
DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( bitand_equal, &= );
DECLARE_ARRAY_BIN_OPS_ON_CONST_VAL( bitor_equal, |= );

template <typename T, unsigned int DIM>
array<T, DIM>& array<T, DIM>::operator>>=( const int shift )
{				
	using binary_assign_ops::right_shift_equal;
	iterator_on_one_element<int> iter_wrap( shift );		
	array_binary_apply<right_shift_equal>::go( begin(), end(), iter_wrap );		
	return *this;
}

template <typename T, unsigned int DIM>
array<T, DIM>& array<T, DIM>::operator<<=( const int shift )
{				
	using binary_assign_ops::left_shift_equal;
	iterator_on_one_element<int> iter_wrap( shift );		
	array_binary_apply<left_shift_equal>::go( begin(), end(), iter_wrap );		
	return *this;
}

//////////////////////////////////////////////////////
//	operator overloading -- binary assign
//	-- on a array
//////////////////////////////////////////////////////

#define DECLARE_ARRAY_BIN_OPS_ON_ARRAY( name, op )			\
template <typename T, unsigned int DIM >				\
array<T, DIM>& array<T, DIM>::operator op( const array<T, DIM>& vec ) 	\
{									\
	using binary_assign_ops::name;					\
	/*array_binary_apply<name>::go( begin(), end(), vec.begin() );*/	\
	array_helper::array_evaluator<name>::evaluate( *this, vec.begin() );	\
	return *this;							\
}

DECLARE_ARRAY_BIN_OPS_ON_ARRAY( plus_equal, += );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( minus_equal, -= );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( multiply_equal, *= );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( divide_equal, /= );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( mod_equal, %= );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( xor_equal, ^= );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( bitand_equal, &= );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( bitor_equal, |= );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( right_shift_equal, >>= );
DECLARE_ARRAY_BIN_OPS_ON_ARRAY( left_shift_equal, <<= );
	
//////////////////////////////////////////////////////
//	operator overloading -- binary logical
//////////////////////////////////////////////////////

#define DECLARE_ARRAY_BIN_LOGIC_OPS_ON_ARRAY( name, op )	\
template <typename T, unsigned int DIM >			\
bool array<T, DIM>::operator op( const array<T, DIM>& vec )	\
{								\
	using binary_ops::name;			\
	return	array_binary_logical_apply<name>		\
		::go( begin(), end(), vec.begin() );		\
}

DECLARE_ARRAY_BIN_LOGIC_OPS_ON_ARRAY( equal_op, == );
DECLARE_ARRAY_BIN_LOGIC_OPS_ON_ARRAY( not_equal_op, != );
/*
DECLARE_ARRAY_BIN_LOGIC_OPS_ON_ARRAY( greater_op, > );
DECLARE_ARRAY_BIN_LOGIC_OPS_ON_ARRAY( less_op, < );
DECLARE_ARRAY_BIN_LOGIC_OPS_ON_ARRAY( greater_equal_op, >= );
DECLARE_ARRAY_BIN_LOGIC_OPS_ON_ARRAY( less_equal_op, <= );
*/

//////////////////////////////////////////////////////
//	operator overloading -- binary assign
//	-- on a expression
//////////////////////////////////////////////////////

#define DECLARE_ARRAY_BIN_OPS_ON_EXPR( name, op )			\
template <typename T, unsigned int DIM >				\
template <class expr_T>									\
array<T, DIM>& array<T, DIM>::operator op( expression<expr_T>&& expr ) 	\
{									\
	using binary_assign_ops::name;					\
	array_helper::array_evaluator<name>::evaluate( *this, expr );	\
	return *this;							\
}

DECLARE_ARRAY_BIN_OPS_ON_EXPR( equal, = );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( plus_equal, += );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( minus_equal, -= );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( multiply_equal, *= );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( divide_equal, /= );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( mod_equal, %= );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( xor_equal, ^= );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( bitand_equal, &= );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( bitor_equal, |= );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( right_shift_equal, >>= );
DECLARE_ARRAY_BIN_OPS_ON_EXPR( left_shift_equal, <<= );

#endif //!< ARRAY_HH